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ABSTRACT 

Motivated by recent observations which detect an outer boundary for starless cores, 
. and evidence for time-dependent mass accretion in the Class and Class I protostellar 

' phases, we reexamine the case of spherical isothermal collapse in the case of a finite 

\ mass reservoir. The presence of a core boundary results in the generation of an inward 

c^d ■ propagating rarefaction wave. This steepens the gas density profile from r~ 2 to r~ 3 or 

steeper. After a protostar forms, the mass accretion rate M evolves through three dis- 
tinct phases: (1) an early phase of decline in M, which is a non-sclf-similar effect due to 
spatially nonuniform infall in the prestellar phase; (2) for large cores, an intermediate 
phase of near-constant M from the infall of the outer part of the self- similar density 
profile; (3) a late phase of rapid decline in M when accretion occurs from the region 
affected by the inward propagating rarefaction wave. Our model clouds of small to 
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intermediate size make a direct transition from phase (1) to phase (3) above. Both the 
\ first and second phase are characterized by a temporally increasing bolometric lumi- 

nosity ibolj while Lboi is decreasing in the third (final) phase. We identify the period 
of temporally increasing Lboi with the Class phase, and the later period of terminal 
accretion and decreasing Lboi with the Class I phase. The peak in Lboi corresponds 
to the evolutionary time when 50% ± 10% of the cloud mass has been accreted by the 
pro tostar . This is in agreement with the classification scheme proposed by Andre et 
al. (^)93). We show how our results can be used to explain tracks of envelope mass 
M anv versus Lboi for protostars in Taurus and Ophiuchus. We also develop an analytic 
formalism which reproduces the protostellar accretion rate. 

Key words: hydrodynamics - ISM: clouds - stars: formation. 
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1 INTRODUCTION 

Recent submillimeter and mid-infrared observations suggest 
that prestellar cores within a larger molecular cloud are 
characterized by a non-uniform radial gas density distribu- 
tion (Ward-Thompson et al. Il999l Bacmann et al. l200Cft . 
Specifically, a flat density profile in the central region of 
size J?fl a t is enclosed within a region of approximately r _1 
column density profile (and by implication an r~ 2 density 
profile) of extent i? m id- Beyond this, a region of steeper den- 
sity (p oc r -3 or greater) is sometimes detected. Finally, 
at a distance 7? e d g e, the column density N seems to merge 
into a background, and fluctuate about a mean value that 
is typical for the ambient molecular cloud. The first two re- 
gions, of extent i?fl at and i? m id, respectively, are consistent 
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with models of unbounded isothermal equilibria or isother- 
mal self-simi lar gr avitation al coll apse (e.g. Chandrasekhar 
Il939l : Larson Il969l : Penston Il969lh In either case, the effect 
of an outer boundary is considered to be infinitely far away 
(i.e. ii m id — > 00 in our terminology). In numerical simula- 
tions of gravitational collapse in which there is a qualitative 
change in the physics beyond some radius (e.g. a transi- 
tion from magnetically supercrit ical t o subcritical mass-to- 
flux ratio: Ciolek & Mouschovias ll99'3 : Basu & Mouschovias 
Il994h . the development of a very steep outer density profile is 
also seen. Finally, larger scale simulations of core formation 
in clou ds wit h uniform background column density (Basu & 
Ciolek 200J) show an eventual merger into a near-uniform 
background column density, demonstrating the existence of 
^edge- The implication of an outer density profile steeper 
than r -3 is that there is a finite reservoir of mass to build 
the star(s), assuming that the gas beyond i? c dg C is governed 
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by the dynamics and gravity of the parent cloud, and thus 
does not accrete on to the star(s) formed within the core. 

An important constraint of the observations are the 
actual sizes of the cores. For example, in the clustered 
star formation regions such as p Ophiuchi protocluster, 
^edge < 5000 AU, and Rcd S c/Rfia.t < 5, while in the more 
extended cores in Taurus, 5000 AU < R cdge < 20 000 AU, 
and 5 < R B d S c/Ra^t J; 10 ( see A ndre et al. Il999l : Andre, 
Ward- Thompson, & Barsony 2000). Clearly, only the latter 
case may approach self-similar conditions. 

Once a central hydrostatic stellar core has formed, the 
mass accretion rate is expe cted t o be co nstant in isother- 
mal similari ty solu tions (Shu ll977t Hunter Il977l Whitworth 
& Summers Il985l) . For example, for the collapse from rest 
of a singular isothermal sphere (SIS) with density profile 
Psis = Cs / (2nGr 2 ), where c s is the isothermal sound speed, 
Shu 1119771) has shown that the mass accretion rate (M) is 
constant and equal to 0.975 c 3 /G. However, two effects can 
work against a constant M in more realistic scenarios of 
isothermal collapse: (1) inward speeds in the prestellar phase 
are not spatially uniform as in the similarity solutions, and 
tend to increase inward, meaning that inner mass shells fall 
in with a greater accretion rate; (2) the effect of a finite mass 
reservoir will ultimately reduce accretion. The first effect has 
been clearly documented in a series of p apers (e.g. Hunter 



119771: Foster & Cheva lier Il993t Tomisaka 11999 : Bas u 11997 
Ciolek & Konigl ll99st Ogino, Tomisaka, & Nakamura[l99j. 

It is always present since the outer boundary condition for 
collapse is distinct from the inner limit of self-similar super- 
sonic infall found in the Larson-Penston solution. Rather, 
the outer boundary condition must represent the ambient 
conditions of a molecular cloud, which d o not correspond 
to large-scale infall (Zuckerman & Evans Il974l) . Addition- 
ally, the finite mass reservoir and steeper than r -3 profile 
as a source of the declining accretion rate has b een stu died 
analytically by Henriksen, Andr e, fc B ontemps ll997l) and 
Whitworth & Ward- Thompson 1200 ll) . although they did 
not account for the physical origin of such a steep density 
slope. 

Indeed, a study of outflow acti vity f rom young stellar 
objects (YSO's) by Bontemps et al. dl99fJ: hereafter BATC) 
suggests that M declines significantly with time during the 
accretion phase of protostellar evolution. Specifically, BATC 
have shown that if the CO outflow rate is proportional to 
M, then Class objects (young protostars at the begin- 
ning of the main accretion phase) have an M that is factor 
of 10 greater (on average) than that of the more evolved 
Class I objects. In this paper, we investigate in detail how 
the assumption of constant mass and volume of a gravita- 
tionally contracting core can affect the mass accretion rate 
and other observable properties after the formation of the 
central hydrostatic stellar core. A very important question 
is: which of the two effects mentioned above - a gradient of 
infall speed in the prestellar phase, or a finite mass reservoir 
and associated steep outer density slope - is more relevant 
to explaining the observations of BATC? The evolutionary 
tracks of envelope mass M env versus bolometric luminosity 
Lboi are another impor tant diagnostic of protostellar evolu- 
tion (Andre et al. l2000f) . BATC have fit the data using a toy 
model in which M decreases with time in exact proportion 
to the remaining envelope mass M env , i.e. M = M cnv /r, 
where r is a characteristic time. 



We seek to explain the observed YSO evolutionary 
tracks using a physical (albeit highly simplified) model. We 
perform high resolution one-dimensional spherical isother- 
mal simulations. The initial peak and decline in the mass 
accretion rate is modeled through numerical simulations and 
a simplified semi-analytic approach. A second late-time de- 
cline in M due to a gas rarefaction wave propagating inward 
from the outer edge of a contracting core, is also studied 
in detail. Comparisons are made with the observationally 
inferred decrease of mass accretion rate (BATC), and evo- 
lutionary tracks of M env ve rsus bolometric luminosity Lboi 
(from Motte & Andre 1200 ll) . 

Numerical simulations of spherical collapse of isother- 
mal cloud cores are described in §|5| The comparison of the 
model with observations is given in § |3| Our main conclu- 
sions are summarized in §0] An analytical approach for the 
determination of the mass accretion rate is presented in the 
Appendix. 



2 ISOTHERMAL COLLAPSE 
2.1 Model Assumptions 

We consider the gravitational collapse of spherical isother- 
mal (temperature T = 10 K) clouds composed of molecu- 
lar hydrogen with a 10% admixture of atomic helium. The 
models actually represent cloud cores which are embedded 
within a larger molecular cloud. The evolution is calculated 
by solving the hydrodynamic equations in spherical coordi- 
nates: 
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where p is the density, v r is the radial velocity, M is the 
enclosed mass, e is the internal energy density and p — (7 — 
l)e is the gas pressure. The ratio of specific heats is equal to 
7 = 1.001 for the gas number density n ^ 10 11 cm~ 3 , which 
implies isothermality (the value of 7 is not exactly unity in 
our implementation in order to avoid a division by zero) . We 
define the gas number density n — p/m, where m — 2.33 ttih 
is the mean molecular mass. When the gas number density 
in the collapsing core exceeds 10 11 cm" 3 , we form the central 
hydrostatic stellar core by imposing an adiabatic index 7 = 
5/3. This simplified treatment of the transition to an opaque 
protostar misses the details of the physics on small scales. 
Specifically, a proper treatment of the accretion shock and 
radiative transfer effects is required to accurately predict 
the p roperties of the stellar core (see Winkler & Newman 
1980 for a detailed treatment and review of work in this 
area). However, our method should be adequate to study the 
protostellar accretion rat e, and has been used succ essfull y by 
e.g. Foster & Chevalier dl993ft and Ogino et al. Il999l) for 
this purpose. We use the method of finite-differences, with 
the time-explicit, operator split solution procedure used in 
the ZEUS-ID numerical hydrod ynam ics code; it is described 
in detail by Stone & Norman 119921) . We have introduced 
the momentum density corre ction factor, as advocated by 
Monchmeyer & Muller (1989), to avoid the development of 
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an anom alous density spike at the origin (see Vorobyov & 
Tarafdar 1999 for details). The numerical grid has 700 points 
which are initially uniformly spaced, but then move with the 
gas until the central stellar core is formed. This provides an 
adequate resolution throughout the simulations. 

We impose boundary conditions such that the gravita- 
tionally bound cloud core has a constant mass and volume. 
The assumption of a constant mass appears to be obser- 
vationally justified by the sharp outer density profiles de- 
scribed in § Physically, this assumption may be justified 
if the core decouples from the rest of a comparatively static, 
diffuse cloud due to a shorter dynamical timescale in the 
gravitationally contracting central condensation than in the 
external region. A specific example of this, due to enhanced 
magnetic support in the outer envelope, is found in the mod- 
els of ambipolar-diffu sion i nduced core formation (see, e.g. 
Basu & Mouschovias Il995ft . The assumption of a constant 
volume is mainly an assumption of a constant radius of grav- 
itational influence of a cloud core within a larger parent dif- 
fuse cloud. 

The radial gas density distribution of a self-gravitating 
cloud with finite central den sity th at is in hydrostatic equi- 
librium (e.g. Chandrasekhar 1939) can be conveniently ap- 
proximated by a modified isothermal sphere, with gas den- 
sity 



(Binney & Tremaine Il987ft . where p c is the central density 
and r c is a radial scale length. We choose a value r c = 
1.1 Cs/\/irGp c , so that the inner profile is close to that of a 
Bonnor-Ebert sphere, r c is comparable to the Jeans length, 
and the asymptotic density profile is 2.2 times the equilib- 
rium singular isothermal sphere value psis = <Z / '(2-kGv 2 ). 
The latter is justified on the grounds that core formation 
should occur in a somewhat non-equilibrium manner (an 
extreme case is the Larson-Penston flow, in which case the 
asymptotic density profile is as high as 4.4 psis), and also by 
observations of protostellar envelope density profiles that are 
often overdense compared to psis (Andre, Motte, & Belloche 
12001ft . 

For small radii (r ^ r c ), the initial density is very close 
to the equilibrium solution for an isothermal sphere with a 
finite central density. However, at large radii it is twice the 
value of the equilibrium isothermal sphere, which converges 
to psis. Hence, our initial con dition s resemble th ose of other 
workers (Foster & Chevalier Il993t Ogino et al. Il999i) who 
start with Bonnor-Ebert spheres and add a positive den- 
sity perturbation in order to initiate evolution. Use of the 
modified isothermal sphere simplifies the analysis a little bit 
since there is a clear transition from flat central region to a 
power-law outer profile. The choice of central density p c and 
outer radius r ou t determines the cloud mass. We study many 
different cloud masses - two models are presented in this sec- 
tion and other models are used to fit observational tracks in 
§ El We also add a (small) positive density perturbation of 
a factor of 1.1 (i.e. the initial gas density distribution is in- 
creased by a factor of 1.1) to drive the cloud (especially the 
inner region which is otherwise near-equilibrium) into grav- 
itational collapse. 

TableQshows the parameters for two model clouds pre- 
sented in this section. The adopted central number density 



Table 1. Model parameters 



Model 


n c 


r c 


r ou t 


Pc 
Pout 




T 


11 


5.0 


0.033 


0.16 


24 


5 


10 


12 


5.0 


0.033 


0.5 


324 


24 


10 



All number densities are in units of 10 4 cm 3 , lengths are in pc, 
masses in Mq, and temperatures in K. 

n c = 5 x 10 4 cm -3 is roughly an order of magnitude lower 
than is observed in prestellar cores (Ward-Thompson et al. 
1999). Considering that these cores may be already in the 
process of slow gravitational contraction, our choice of n c is 
justified for the purpose of describing the basic features of 
star formation. In both models, the outer radius r out is cho- 
sen so as to form gravitationally unstable prestellar cores 
with central-to-surface density ratio p c /p ou t ^ 14 (since 
our initial states are similar to Bonnor-Ebert spheres). In 
model II, Pc/pout ~ 24 and by implication r out /r c w 4.7, 
whereas in model 12 p c /p ou t ~ 324 and r out /r c w 18. 
Model 12 thus represents a very extended prestellar core. 
Models II and 12 have masses 5 Mq and 24 Mq respec- 
tively; the T stands for isothermal. 



2.2 Numerical Results 

Fig.^shows the temporal evolution of the radial gas density 
profiles (the upper panel) and velocity profiles (the lower 
panel) during the runaway collapse phase (before the for- 
mation of the central hydrostatic stellar core) in model II. 
The density and velocity profiles are numbered according to 
evolutionary sequence, starting from the initial distributions 
(profile 1; note that the cloud core is initially at rest) and 
ending with those obtained when the central number den- 
sity has almost reached 10 10 cm -3 (profile 5). The dashed 
lines in the upper panel of Fig. show the power-law in- 
dex din p/dlnr of the gas distribution for profiles 1, 2, and 
3. By the time that a relatively mild center-to-boundary 
density contrast ~ 150 is established (profile 2), the radial 
density profile st arts re sembling those observed in Taurus by 
Bacmann et al. i2000ft : it is flat in the central region, then 
gradually changes to an r~ 2 profile, and falls off as r~ 3 or 
steeper in the envelope at r > 0.08 pc. The sharp change in 
slope of the density profile (e.g. at r ~ 0.08 pc in profile 2 
of Fig. is due to an inwardly-propagating gas rarefaction 
wave caused by a finite reservoir of mass. The self-similar re- 
gion with r~ 2 density profile is of the Larson-Penston type, 
with density somewhat greater than the equilibrium singu- 
lar isothermal sphere value (c 2 /2-KGr 2 ). The velocity profiles 
in Fig. also show a distinct break at the instantaneous lo- 
cation of the rarefaction wave. Furthermore, the peak infall 
speed is clearly supersonic (since c s = 0.19 km s _1 ) by the 
time profile 4 is established, again consistent with Larson- 
Penston type flow in the inner region. 

Fig. |2K shows the temporal evolution of the accretion 
rate at a radial distance of 600 AU from the center in 
model II. The evolution is characterized by a slow initial 

1 We note that the accretion rate is not expected to vary signifi- 
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Figure 1. The radial gas density (the upper panel) and velocity 
(the lower panel) profiles obtained in model II before the central 
hydrostatic stellar core is formed. The number I corresponds to 
the initial profiles (note that initially the cloud core is at rest) 
and the number 5 labels profiles when the central gas number 
density has almost reached 10 10 cm -3 . The dashed lines in the 
upper panel show the power-law index dlnp/dlnr of the gas dis- 
tributions 1, 2, and 3. Profiles 2, 3, 4, and 5 are reached at times 
0.309 Myr, 0.378 Myr, 0.392 Myr, and 0.394 Myr, respectively. 
For reference, the dotted line is the density profile of a singular 
isothermal sphere. 



gravitational contraction and then a very rapid increase un- 
til about 0.4 Myr. Subsequently, a central hydrostatic stellar 
core forms and the mass accretion rate reaches its maximum 
value of 2.8 x 10~ 5 M Q yr _1 (or 17.4 4/G). After stellar 
core formation, the evolution of the mass accretion rate has 
possibly three distinct phases, of which two are on display 
m Fig. Hi. The early phase, plotted with the dashed line 
in Fig. is characterized by accretion of material that 
has not yet been affected by the rarefaction wave propa- 
gating inward from the outer boundary. The accretion rate 
is declining, even though the density profile near the stel- 
lar core was nearly self-similar at the moment of its forma- 
tion. This decline is due to the gradient of infall velocity 
in the inner regions, an effect not predicted in the similar- 
ity solutions. However, if there is a large outer region with 
mass shells that are falling in at significantly subsonic speeds 
when the central stellar core forms (see discussion of Fig. 



cantly in the range 0.1 AU < r < 1000 AU according to Masunaga 
& Inutsuka <2000h. 





\ (65.0) 

\ 










\ 

\ 

\ 








- 


. (40.0) 






- 


[ , 
J 


\ 


(5.0) 






; 














> 

\ (93.0) 

\ 


1 








\ 

vl 

\ 

\ 

\ 

\ 




(37.0) 




u 








(50) 























0.3 0.5 



Time (Myr) 



Figure 2. The temporal evolution of the mass accretion rate at 
the radial distance r = 600 AU from the center of a cloud core 
obtained in a) model II and b) model 12. The model cloud II has 
mass 5 Mq and the model cloud 12 has mass 24 Mq. The solid 
lines show M during the runaway collapse phase, prior to the for- 
mation of a central stellar core. The dashed and dotted lines plot 
M after stellar core formation; the dashed lines show M before 
the gas affected by the inwardly propagating rarefaction wave has 
reached r = 600 AU, whereas the dotted lines show M after this 
gas has reached r = 600 AU. The numbers in parentheses reflect 
the percentage of cloud mass remaining in the envelope at the 
given times. 



below), the accretion rate will eventually stabilize to a 
constan t valu e that is consistent with the standard theory 
of Shu lll977l) . In that picture, progressively higher shells 
of gas lose their partial pressure support and start falling 
from rest on to the central stellar core almost in a free-fall 
manner. This would be the intermediate phase of accretion. 
However, the late phase of very rapid decline of the accre- 
tion rate starts at roughly 0.46 Myr (before the intermediate 
phase can be established in the 5 Mq cloud), when gas af- 
fected by the inwardly propagating rarefaction wave reaches 
the inner 600 AU. This results in a sharp drop of M as shown 
in Fig. [5^ by the dotted line. 

The existence of the (in principle) three distinct phases 
of mass accretion is clearly seen Fig. |!J>, where M of the 
more extended cloud (r out /r c w 18) is plotted (hereafter, 
model 12). The outer boundary is now at r out = 0.5 pc and 
it takes a time > 1 Myr for the influence of the rarefac- 
tion wave to reach the inner 600 AU. As a result, the mass 
accretion rate has time to stabilize at a constant value of 
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1.34 x 1CT 5 M Q yr _1 (the dashed line in Fig. Us), before it 
sharply drops a t later times (the dotted line in Fig. |5p) . Ac- 
cording to Shu the collapse from rest of a power-law 
profile that has a density equal to twice psis yields a mass 
accretion rate 5.58 c 3 /G = 8.86 x 1CT 6 M Q yr" 1 . Our stable 
intermediate accretion rate is roughly consistent with this 
prediction since the density in the power-law tail is actually 
somewhat greater than twice psis- It is equal to 2.42 psis 
in the initial state, and grows to greater overdensities in 
the innermost regions. However, the bulk of the matter, 
which is in the outer tail, has density within 2.5psis- Fur- 
ther experiments with our numerical simulations show that 
the intermediate phase of constant accretion rate is observed 
only in rather exten ded p restellar cores with r out /r c > 15. 
Foster & Chevalier |12S3) found an even stronger criterion 
r ut/r c > 20. Since more extended cores tend to be more 
massive as well, we may expect to observe the intermediate 
phase more frequently in the collapse of massive cores. 



2.3 Effect of Boundary Condition 

Our standard simulation does not contain an external 
medium explicitly. In order to explore the effect of such a 
medium, we ran additional simulations in which the cloud 
core is surrounded by a spherical shell of diffuse (i.e. non- 
gravitating) gas of constant temperature and density. The 
outermost layer of the cloud core and the external gas are 
initially in pressure balance. We found that the value of M 
in the late accretion phase may depend on the assumed val- 
ues of the external density and temperature. For instance, 
if the gravitating core is nested within a larger diffuse non- 
gravitating cloud of T = 10 K and p = p ou t, the accretion 
rate increases slightly as compared to that shown in Fig.|5|by 
the dotted line. A warmer external non-gravitating environ- 
ment of T — 200 K and p = p out /20 shortens the duration of 
the late accretion phase shown in Fig. [3] by the dotted line. 
This phase may be virtually absent if the sound speed of the 
external diffuse medium is considerably higher (by a factor 
~ 1000) than that of the gravitationally bound core. This 
essentially corresponds to a constant outer pressure bound- 
ary condition (see Foster & Chevalier [1993). However, such 
a high sound speed contrast is not expected for star forma- 
tion taking place in a dense (n ~ 10 4 cm" 3 ) environment 



like p Ophiuchi (Johnstone et al. 12000 

We believe that the constant volume boundary condi- 
tion, and resulting inward propagating rarefaction wave, are 
best at reproducing the steep outer density profiles and the 
low (residual) mass accretion rate necessary to explain the 
Class I phase of protostellar accretion. 



2.4 Semi-analytic Model 

Finally, we compute M of a pressure-free cloud using the 
analytical approach developed in the Appendix. This ap- 
proach allows for the determination of M for a cloud with 
given initial radial density p(ro) and velocity vo(ro) pro- 
files, if the subsequent collapse is pressure-free. We find that 
the success or failure of the analytical approach to describe 
the mass accretion rate of the isothermal cloud depends on 
the adopted p(ro) and vo(ro) profiles. For instance, if p(ro) 
is determined by profile 2 (the upper panel of Fig. and 
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Figure 3. The temporal evolution of the mass accretion rate. The 
solid line shows the results of isothermal numerical simulations 
(model II). The dashed line shows M obtained in a pressure-free 
approximation if such collapse begins from rest with the relatively 
mildly concentrated density profile 2 in Fig. The dotted line 
shows the result for pressure-free collapse if a non-zero initial 
velocity (that of profile 2 in the bottom panel of Fig. is also 
used. The agreement of the pressure- free model and full numerical 
simulation are quite good in the latter case. 



vo{ro) = 0, the pressure- free mass accretion rate shown in 
Fig. |3] by the dashed line reproduces only very roughly the 
main features of the isothermal accretion rate (the solid line 
in Fig. El . However, if we take into account the non-zero and 
non-uniform velocity profile vo(ro) plotted in the lower panel 
of Fig. [T] (profile 2), then the pressure-free M shown by the 
dotted line in Fig. ^reproduces that of the isothermal cloud 
much better. This example demonstrates the importance of 
the velocity field prior to stellar core formation in deter- 
mining the accretion rates after its formation. The success 
of our analytical pressure-free approach also shows that the 
collapse of the isothermal cloud can be regarded as essen- 
tially pressure- free from the time of a relatively mild central 
concentration p c /pout > 150, when the central number den- 
sity ~ 2 x 10 5 cm -3 . 



3 ASTROPHYSICAL IMPLICATIONS 

Class objects represent a very e arly phase of protostellar 
evolution (see Andre et al. 1200(f ) . as evidenced by a rela- 
tively high ratio of submillimeter luminosity to bolometric 
luminosity: L BUOlnnl / L no \ > 0.5%. Class objects also drive 
powerful collimated CO outflows. A study of outflow ac- 
tivity in low-mass YSO's by BATC suggests that the CO 
momentum flux F co declines significantly during protostel- 
lar evolution. Specifically, F co decreases on average by more 
than an order of magnitude from Class to Class I objects. 
This tendency is illustrated in Fig. |1J where we plot F co ver- 
sus M onv for 41 sources listed in BATC. We relate F co to M 
by 



F co = .font X (M w /M)Vw x M, 



(•») 



where / cnt is the entrainment efficiency that relates F co to 
the momentum flux M W VW of the wind. Based on theoretical 
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Figure 4. CO momentum flux F co versus envelope mass M en v 
for 41 sources of Bontemps et al. (1996). The Class and Class I 
objects are plotted with the open and filled circles, respectively. 
The model F co — M env tracks of three prestellar clouds of M c \ = 
0.2 Mq, 0.7 Mq, and 2 Mq are shown by the dotted, dashed, 
and solid lines, respectively. 



Table 2. Model parameters for Ophiuchus 



M cl 


n c 


r-out 


Pc/Pout 


r ou t /r c 


a 


0.17 


1 x 10 7 


1600 


15.0 


3.7 


2.0 


0.23 


5 x 10 6 


1900 


14.1 


3.6 


2.0 


0.55 


2 x 10 6 


4000 


18.0 


i.l 


2.1 


0.9 


2 x 10 6 


4000 


18.0 


i.l 


4.0 



All number densities are in cm 3 , lengths in AU, and masses in 

Mq. 



Table 3. Model parameters for Taurus 





n c 


r ou t 


Pc/Pout 


r ou t/r c 


a 


0.46 


1.5 x 10 6 


5000 


18 


4.1 


1.9 


0.65 


1 x 10 6 


6000 


26 


5.0 


1.8 


1.0 


1 x 10 6 


10000 


71 


8,1 


1.5 


1.5 


8 x 10 5 


12000 


73 


8.5 


1.8 



All number densities are in cm 3 , lengths in AU, and masses in 
Mq. 



models in the literature, BATC suggested that the factor 
/cut and the outflow driving engine efficiency (M w /M)Vw 
do not vary significantly during protostellar evolution. This 
implies that the observed decline of _F co reflects a corre- 
sponding decrease in M from the Class to the Class I 
stage. Following BATC, we take / ont = 1, M w /M = 0.1, 
and Vw = 150 km s _1 and use equation JKJ to compute F co 
from our model's known mass accretion rate M (see Fig.|5|. 

Since the sample of Class and Class I objects listed in 
BATC includes sources from both the Ophiuchus and Tau- 
rus star forming regions, we develop model clouds which 
take into account the seemingly different initial conditions 
of star formation in these regions. As mentioned in § Q 
the two most prominent differences between these two re- 
gions are: (1) The cores in Ophiuchus have outer radii 
(r ou t < 5000 AU) which are smaller than in T aurus , where 
5000 AU < r ut < 20000 AU (Andre et al. Il999l: Andre 
et al. 12(10(1; ; (2) The radial column density profiles of the 
protostellar envelopes of Class objects in Ophiuchus are 
at least 2-3 times denser than a SIS at T — 10 K, whereas in 
Taurus the protostellar envelopes are overdense co mpare d 
to the SIS by a smaller factor < 2 (Andre et al. 1200 ll) . 
This implies that radial column density profiles of prestellar 
cores in Ophiuchus and Taurus may follow the same ten- 
dency. We develop a set of Ophiuchus model cores which 
have r ou t < 5000 AU, and a set of Taurus model cores which 
have 5000 AU < r out < 20000 AU. Furthermore, the factor 
a (by which our model density profiles are asymptotically 
overdense compared to psis) is taken to be > 2.0 for Ophi- 
uchus and < 2.0 for Taurus. Clearly, there is no unique set 
of model cloud parameters that would be exclusively con- 
sistent with the observational data, given the measurement 
uncertainties. We have chosen a set of core central densities 
p c , radii r out , and overdensity factors a so as to reasonably 
reproduce the observed properties of the cores in the two 
regions. The parameters of the model density distributions 
for Ophiuchus and Taurus are listed in Tableland Tabled 
respectively. 



We also note that we have ensured that the cores satisfy 
the gravitational instability criterion r out /r c > 3.6, which is 
similar to that for Bonnor-Ebert spheres. The Ophiuchus 
model cores are clustered near this limiting value of r out /r c , 
but the Taurus model cores are allowed to be somewhat 
more extended, again in keeping with observed properties. 
We also note that the masses of prestellar cores with the 
radial density profile given by equation Q scale as 1/pS?' 5 , 
if the ratio r out /r c is fixed. 

The sample of 41 sources in BATC contains Class and 
Class I objects from both Ophiuchus and Taurus. Hence, 
in Fig. 0] we take three representative prestellar clouds of 
M c i = 0.23 Mq (Ophiuchus), 0.65 Mq (Taurus), and 2.0 Mq 
(Taurus), for which the F co — M cnv tracks are shown by 
the dotted, dashed, and solid lines, respectively. Both the 
data and model tracks show a near-linear correlation be- 
tween F co (oc M) and M onv . A slightly better fit of the model 
tracks to the data can be obtained by adjusting one or more 
of the estimated parameters / cn t, M w /M, and K» by factors 
of order unity. 

Based on the near- linear correlation of F co and M env , 
BATC developed a toy model in which M decreases with 
time in exact proportion to the remaining envelope mass 
M env , i-e. M = Menv/r, where r is a characteristic time. 
Furthermore, if one assumes that the bolometric luminos- 
ity derives entirely from the accretion on to the hydrostatic 
stellar core, i.e., L DO i = GM C M / R c , where M c and 7i c are 
the mass and radius of the stellar core, respectively, then 
the bolometric luminosity reaches a maximum value when 
half of the initial prestellar mass has been accreted by the 
protostar and the other half remains in the envelope. The 
evolut ionary time when M c = M env was defined by Andre 
et al. ( 1993?) as the conceptual border between the Class 
and Class I evolutionary stages. 

The solid and dashed lines in Fig. [HJi and Fig. |H|d show 
Lboi and F co obtained in model II and model 12, respec- 
tively. Since we do not follow the evolution of a proto- 
star to the formation of the second (atomic) hydrostatic 
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Time (Myr) 



-bol 



(L, 



class class I 





Time (Myr) 

Figure 5. Temporal evolution of the bolometric luminosity Lbol 
and CO momentum flux F co . a) The solid and dashed lines show 
Lbol an d F co obtained in model II, respectively, b) The solid and 
dashed lines show Lbol and _F C o obtained in model 12, respectively. 
Note that Lbol is still increasing during the early phase of accre- 
tion rate decline but only declines later due to the more severe 
accretion rate decline caused by the inward propagating rarefac- 
tion wave. The vertical dash-dotted line is the temporal dividing 
line between the Class and Class I phases for each model; the 
numbers below give the mass of the central stellar core as a per- 
centage of the total cloud mass. Crosses indicate the time when 
50% of the initial cloud mass has been accreted by the protostar. 
Note the use of a logarithmic scale for time, so that the Class I 
phase is still longer than the Class phase for model 12. 



1.5 M 



1.0 Mg _ _ 
0.65 M . _ • • 

0.46 M 




(a) Taurus 



0.9 Ma, 



0.55 M, 

0.23 M @ 
0.17 M .. - . _T 




(b) Ophiuchus 



L b0 i (Le) 

Figure 6. Envelope mass M cn v versus bolometric luminosity Lboi 
for 37 protostellar objects taken from Mottc & Andre (2001). Di- 
agrams are shown for a) Taurus, and b) Ophiuchus. The Class 
and Class I objects are plotted with the open and filled circles, 
respectively. The triangles represent the observed peculiar Class I 
sources. The model M cnv — Lbol tracks of eight prcstcllar clouds 
with masses M cl = 1.5 Mq, 1.0 Mq , 0.65 Mq, and 0.46 Mq (Tau- 
rus) and 0.9 Mq, 0.55 Mq, 0.23 Mq, and 0.17 Mq (Ophiuchus) 
are shown by the solid, dashed, dotted, and dotted-dashed lines, 
respectively. Crosses indicate the time when 50% of the cloud 
mass has been accreted by the protostar. 



core, we take R c = 3 Rq and let Lboi = GM C M / R c . 
The radius R c depends on th e acc retion rate and stel- 
lar mass (see Fig. 7 of Stahler Il98ctl and may vary from 
R c « 1.5 Rq for small stellar cores M c ~ 0.2 Mq and low 
accretion rates M~2x 10" 6 Mq yr" 1 to Rc « 5.0 Rq for 
large stellar cores M c ~ 1.0 Mq and high accretion rates 
M ~ 1 x 10 _J Mq yr -1 . However, this variation consti- 
tutes roughly a factor of 2 change in the adopted average 
value of R c — 3 Rq. Indeed, we performed numerical sim- 
ulations with a varying R c (assuming a normal deuterium 
abundance) and found that it has only a minor qualitative 
influence on our main results. The stellar core mass M c is 
computed by summing up the masses of the central hydro- 
static spherical layers in our numerical simulations. An ob- 
vious difference in the temporal evolution of F co and Lboi is 
seen in Fig. [S] The temporal evolution of F co after the central 
hydrostatic core formation at t w 0.4 Myr goes through the 
same phases as shown for M in Fig. The temporal evolu- 
tion of Lboi shows two distinct phases: it increases during the 
early phase (unlike F co ) and starts decreasing only when gas 



affected by the inward propagating rarefaction wave reaches 
the central hydrostatic core. Thus, in our model, only the 
rarefaction wave acts to reduce Lboi during the accretion 
phase of protostellar evolution. This is a physical explana- 
tion for the peak in Lboi that also occurs in the toy model 
of BATC. In that model, the bolometric luminosity reaches 
a maximum value when exactly half of the initial prestellar 
mass has been accreted by the protostar. In our simulations, 
the peak in Lboi corresponds to the evolutionary time when 
50% ± 10% of the matter is in the protostar (higher devia- 
tions up to +15% are found in very massive and extended 
prestellar clouds). 

Finally, in Fig. |S| we show the M onv — Lboi evolution- 
ary tracks. We use eight representative prestellar cloud core 
masses as listed in Table H and Table Fig. |SJi shows the 
overlaid data for YSO's in Taurus, while Fig. |S|d has over- 
laid data for Ophiuchus. Th e dat a for both samples are 
taken from Motte & Andre (20QlJ). The open circles rep- 
resent bonafide Class objects, the solid circles represent 
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the bonafide Class I objects, while the triangles represent 
the so-called peculiar Class I objects observed in Taurus. 
We note t hat th e envelope masses given in Table 2 of Motte 
& Andre (2001) and plotted in their Fig. 5 and Fig. 6 are de- 
termined within a 4200 AU radius circle. While this should 
relatively well describe the total envelope masses in Ophi- 
uchus, a substantial (a factor of 3) portion of the envelope 
mass may be missing in the Taurus cores, which have sizes as 
large as 15000-20000 AU. For this reason, we plot in Fig.EJi 
the to tal envelope masses given in Table 4 of Motte & Andre 
(2001) for a set of resolved Taurus cores. 

The loci of maximum Lboi in the M em — Lboi tracks 
roughly separate two phases in the evolution of a protostar: 
a shorter one characterized by accretion of matter from the 
envelope not yet affected by the rarefaction wave (i.e. char- 
acterized by the r~ 2 gas density profile or shallower) and 
a longer one characterized by accretion of matter from the 
rarefied envelope (i.e. characterized by the r -3 profile or 
steeper). The turnover also corresponds to the evolutionary 
time when 50% ± 10% of the matter is in the protostar and a 
corresponding amount remains in the envelope as shown by 
the crosses in Fig. © This is in agreement with the observa- 
tional requirements and toy model of BATC. Given that the 
peak in Lboi is our conceptual dividing line between two dis- 
tinct phases of accretion, we conclude that in Taurus, most 
of the so-called Class I objects would tend to fall into the 
Class category in our scheme. They may indeed be more 
evolved than the already identified Class objects, having 
lower values of M and M cnv , but would not be in a qualita- 
tively distinct phase of evolution (see Motte & Andre EoOll 
for a similar conclusion). In contrast, the so-called peculiar 
Class I objects in Taurus would be proper Class I objects in 
our scheme since they are likely in a phase of declining Lboi- 
In Ophiuchus, the currently identified Class and Class I 
objects do seem to fall on two distinct sides of the peak in 

Lboi- 

Fig.|^]indicates that in extended clouds (as of model 12) 
the phase of increasing Lboi is longer than in compact clouds 
(as of model II). This may explain why this phase is more 
populated in Taurus than in Ophiuchus. In addition, ex- 
tended clouds have a longer phase of accretion from the en- 
velope not yet affected by the rarefaction wave and, as a con- 
sequence, a higher probability of having a quasi-constant ac- 
cretion phase. This is in agreement with the previous sugges- 
tions made by Henriksen et al. il997f) and Andre et al. ( 2000) 
that the accretion history in Taurus is closer to the SIS sce- 
nario than in Ophiuchus. However, we note that none of our 
representative prestellar cores listed in Table [3] and used to 
fit the data for Taurus are extended enough (r ou t/r c > 15) 
to have a distinct phase of constant accretion. 

One problem should be pointed out here. While our 
model Lboi — M cnv tracks in Fig. |S| explain well the mea- 
sured bolometric luminosities in Ophiuchus, they seem to 
overestimate Lboi in Taurus by a factor of 5-10. This is 
the so-called "luminosity proble m" th at was first noticed 
by Kenyon, Calvet, & Hartmann |l993). As a consequence, 
the position near the turnover in Lboi — M cnv tracks for Tau- 
rus is scarcely populated. This implies that while spherical 
collapse models may be appropriate for the determination 
of Lboi in Ophiuchus, they tend to overestimate L bo i in Tau- 
rus. It is possible that a significant magnetic regulation of 
the early stages of star formation in Taurus, as implied by 



e.g. polarization maps (Moneti et al. ll984l) would yield more 
flattened envelopes which result in a lower accretion rate on 
to the central protostar and a s maller bolometric luminos- 
ity. Interestingly, Kenyon et al. 1^93) also concluded that 
envelopes in Taurus should be highly flattened in order to 
explain their spectral energy distribution. Two dimensional 
simulations are required to address this issue. 

Finally, it is worth noting that our Taurus model cores 
are generally more massive than the Ophiuchus model cores. 
This in agreement with observations, and can be justified 
theoretically on the basis of a lower mean column density 
(hence greater Jeans length and Jeans mass of a sheetlike 
configuration) in Taurus compared to regions of more clus- 
tered star formation in e.g. Ophiuchus and Orion. Taken at 
face value, our models then imply that Taurus protostars 
should be more massive in general than Ophiuchus proto- 
stars. While such a conclusion must be tempered by the 
fact that we do not model magnetic support or feedback 
from outflows, there is some evidence that Taurus does have 
a significantly higher mass peak in its initial mass fun ction 
than does the Trapezium cluster in Orion (see Luhman l2004l 
and references within). 



4 CONCLUSIONS 

Our numerical simulations indicate that the assumption of 
a finite mass reservoir of prestellar cores is required to ex- 
plain the observed Class to Class I transition. We start 
our collapse calculations by perturbing a modified isother- 
mal sphere profile (eq. ^]) that is truncated and resembles a 
bounded isothermal equilibrium state. Specifically, we find 
that 

• Starting in the prestellar runaway collapse phase, a 
shortage of matter developing at the outer edge of a core gen- 
erates an inward propagating rarefaction wave that steepens 
the radial gas density profile in the envelope from r -2 to r~ 3 
or even steeper. 

• After a central hydrostatic stellar core has formed, and 
the cloud core has entered the accretion phase, the mass ac- 
cretion rate M on to the central protostar can be divided 
into three possible distinct phases. In the early phase, M de- 
creases due to a gradient of infall speed that developed dur- 
ing the runaway collapse phase (such a gradient is not pre- 
dicted in isothermal similarity solutions). An intermediate 
phase of near-constant M follows if the core is large enough 
to have an extended zone of self-similar density profile with 
relatively low infall speed during the prestellar phase. Fi- 
nally, when accretion occurs from the region affected by the 
inward propagating rarefaction wave, a terminal and rapid 
decline of M occurs. 

• A pressure-free analytic formalism for the mass ac- 
cretion rate can be used to predict the mass accretion rate 
after stellar core formation, given the density and velocity 
profiles in a suitably late part of the runaway collapse phase. 
Our formulas can estimate M at essentially any radial dis- 
tance from the central singularity. This makes it possible to 
obtain M as a function of radial distance at any given time. 
We have demonstrated the importance of the velocity field 
of a collapsing cloud in determining M; our approach suc- 
cessfully estimates the accretion rate if the velocity field is 
taken into account. It demonstrates that the initial decline 
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in M is due to the gradient of infall speed in the prestellar 
phase. 

• From an observational point of view, we can under- 
stand evolutionary Mew — Lboi tracks using core models 
of relatively small mass and size, so that there is not an 
extensive self-similar region, in a greem ent with the profiles 
observed by e.g. Bacmann et al. 1)200011 . This means that in 
the accretion phase, M makes a direct transition from the 
early decline phase to the late decline phase when matter 
is accreted from the region of steep (r~ 3 or steeper) density 
profile that is affected by the inward propagating rarefac- 
tion wave. In the first phase (which we identify as the true 
Class phase), the bolometric luminosity Lboi is increasing 
with time, even though M and the CO momentum flux L co 
are slowly decreasing. In the second phase (which we identify 
as the Class I phase), both Lboi and F co decline with time. 
Hence, our simulations imply that the influence of the rar- 
efaction wave roughly traces the conceptual border between 
the Class and Class I evolutionary stages. Regions of star 
formation with more extended cores, like Taurus, should re- 
veal a larger fraction of protostars in the phase of increasing 
Lboi- Our Fig. [5] reveals that this is indeed the case, if most 
of the so-called Class I objects in Taurus are reclassified as 
Class 0, according to our definition. The so-called peculiar 
Class I objects in Taurus would be bona-fide Cla ss I ob jects 
according to our definition (see Motte & Andre EoQll for a 
similar conclusion on empirical grounds). 

• Luminosities derived entirely from the accretion on 
to the hydrostatic stellar core tend to be larger than the 
measured bolometric luminosities Lboi in Taurus by a fac- 
tor of 5-10, while they seem to better explain the measured 
Lboi in Ophiuchus. This implies that physical conditions in 
Ophiuchus may favour a more spherically symmetric star 
formation scenario. 

Our results should be interpreted in the context of mod- 
els of one-dimensional radial infall. They illuminate phenom- 
ena which are not included in standard self-similar mod- 
els of isothermal spherical collapse, by clarifying the impor- 
tance of boundary (edge) effects in explaining the observed 
Leo — Menv and M cnv — Lboi tracks. Important theoretical 
questions remain to be answered, such as the nature of the 
global dynamics of a cloud which could maintain a finite 
mass reservoir for a core. A transition to a magnetically 
subcritical envelope may provide the physical boundary that 
we approximate in our model. For example, Shu, Li, & Allen 
(2004) have recently calculated the (declining) accretion rate 
from a subcritical envelope on to a protostar, under the as- 
sumption of flux freezing. An alternate or complementary 
mechanism of limiting the available mass reservoir is the 
effect of protostellar outflows. 

Our main observational inference is that a finite mass 
reservoir and the resulting phase of residual accretion is nec- 
essary to understand the Class I phase of protostellar evo- 
lution. Our calculated mass accretion rates really represent 
the infall onto an inner circumstellar disk that would be 
formed due to rotation. Hence, our results are relatable to 
observations if matter is cycled through a circumstellar disk 
and on to a protostar rapidly enough so that the protostel- 
lar accretion is at least proportional to the mass infall rate 
on to the disk. This is likely, since disk masses are not ob- 
served to be greater than protostellar masses, but needs to 
be addressed with a more complete model. In future papers, 



we will investigate the role of non-isothermality (using de- 
tailed cooling rates due to gas and dust), rotation, magnetic 
fields, and non-axisymmetry in determining M and implied 
observable quantities. 
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APPENDIX A: PRESSURE-FREE COLLAPSE 
Al Collapse from rest 

The equation of motion of a pressure- free, self-gravitating 
spherically symmetric cloud is 

dv _ GM (r) 
dt ~ r^~ ' 
where v is the velocity of a thin spherical shell at a radial 
distance r from the center of a cloud, and M(r) is the mass 
inside a sphere of radius r. equation 1AH can be integrated 
to yield the expression for velocity v 



(Al) 



dr 



2GM(r ) 



I _ _L 

r r 



(A2) 



where ro is the initial position of a mass shell at t — 0, and 
M(ro) is the mass inside ro. Here, it is assumed that all 
shells are initially at rest: Vo(ro) = 0. Equation l|A20 can 
be integrat ed by means of the substitution r — ro cos 2 p 
(see Hunter |l962j) to determine the time it takes for a shell 
located initially at ro to move to a smaller radial distance r 
due to the gravitational pull of the mass M(ro). The answer 
is 



i yr/ro + 0.5 sin(2 arccos y/r/ro) 
y/2GM(r )/ti ' 



(A3) 



The velocity v(r, t) at a given radial distance r and time t can 
now be obtained from equations 1A2I and IA3I . The values 
of r and t are sufficient to determine ro (a value > r but ^ 
r ou t, where r out is the radius of a cloud) from equation 1A3I . 
Subsequently, we use the obtained value of ro in equation 
<A2> to obtain v(r,t). 

Provided that the shells do not pass through each other 
(i.e. the mass of a moving shell is conserved, dM(r,t) = 
dM(ro,to)), the gas density of a collapsing cloud is 



p(r,t) = 



Po(ro) rg dr 
r 2 dr 



(A4) 



where po(fo) is the initial gas density at ro. The ratio of 
dro / dr determines how the thickness of a given shell evolves 
with time. The relative thickness dro/dr is determined by 
differentiating r = ro cos 2 (3 with respect to ro, yielding 



dr 
dr Q 



ro 



ro sin 2 arccos 



r \ d/3 
ro ) dr ' 



(A5) 



Next, d/3/dro is determined from an alternate form of equa- 
tion l|A3|l : 



P + 0.5 sin 2/3 = t 



2GM(r ) 



Differentiating with respect to ro yields 



dp 
dr 



G ( dM(r ) 
2M{r )rl I dr 



3M(r ) 



ro 



(A6) 



(A7) 



Now that the density p(r, t) and velocity v(r, i) distribu- 
tions of a collapsing pressure-free sphere are explicitly deter- 
mined, the mass accretion rate at any given radial distance 
r and time t can be found as M(r,t) — iivr 2 p(r,t)v(r,t). 

A2 Collapse with non-zero initial velocity 

In a general case of non-zero initial radial velocity profile 
vo(ro), integration of eqaution (IAH yields 



dr 
dt 



2GM(r ) 



1 

ro 



(A8) 



where uo(ro) is the initial velocity of a shell at ro. Equa- 
tion 1A8II can be reduced to an integrable one by means of 
a substitution r = ro cos 2 P 



(A9) 



sin(2/3) dp = dt 2GM{rp) 
^tan 2 P + a V r o 

where a = ro«o(ro)/2GAf(ro). Another substitution 
sin 2 P — x and integration over x from x = to x — sin 2 P 
finally gives the time t it would take for a shell located ini- 
tially at ro and having a non-zero initial velocity vo to move 
to a smaller radial distance r: 



t = 



(l-a)2GM(r ) 



(l-f + 5) (f 
V r / \r 



-Vd- (1 + 5) tan" 1 V5 



(1+5) tan 



_j y/l- r/r + 8 



(A10) 



where 8 — o/(l — a). 

In the case of a non-zero initial velocity profile, it is 
more complicated to obtain a simple analogue to equations 
IA5I - IA7I and explicitly determine a density distribution 
p(r,t), as done in the previous example. Instead, we obtain 
the mass accretion rate by computing the mass that passes 
the sphere of radius r during time At, i.e. 

A/(r + Ar )-M(r ) 



M(r,t) = 



At 



(All) 



where M(ro) is the mass inside a sphere of radius ro. A time 
interval At is the time that it takes for two adjacent shells 
of radius ro and ro + Aro to move to the radial distance r. 
The value of At can be found by solving equation HAIOII for 
fixed values of ro, ro + Aro, and r. 




Figure Al. The mass accretion rate as a function of time and 
radial distance from the center of a pressure-free cloud that has 
initial radial gas density distribution p = p c /[l+ (r/r c ) 2 ], in which 
pc = 7.5 X 10 4 cm" 3 and r c = 0.033 pc. 




Figure A2. The same as Fig. lAll but for the initial gas density 
distribution p = p c /[l + (f/?*c) 3 ]- 



A3 Applications 

As two examples, we consider two different initial gas density 
profiles and determine the pressure-free mass accretion rate 
M(r,t) as a function of radial distance r and time t. 



A3.1 Modified isothermal sphere 

First, we consider the radial gas density profile of a mod- 
ified isot herma l sphere: p = p c /[l + {r/r c ) 2 ] (Binney & 
Tremaine 0.987) , where p c is the gas density in the center 
of a cloud and r c is the radial scale length. Figure lATl shows 
M(r,t) of a pressure-free cloud with p c — 5.5 x 10 4 cm -3 
and r c = 0.033 pc. The mass accretion rate increases with 
time and appears to approach a constant value at later 
times t > 0.7 Myr. Note that the temporal evolution of 
the mass accretion rate depends on the radial distance r: M 
approaches faster a constant value at smaller r. This behav- 
ior of M(r,t) is independent of the adopted values of p c and 



A3. 2 A steeper profile 

The submillimeter and mid-infrared observatio ns of Ward- 
Thompson et al. (1999) and Bacmann et al. (2000) sug- 
gest that the gas density in the envelope of a starless core 
falls off steeper than r~ 2 . As a second example, we con- 
sider a pressure-free cloud with the initial gas density profile 
p = p c /[l + (r/r c ) 3 ] and plot the corresponding mass accre- 
tion rate M(r,t) in Fig. IA1I The values of p c and r c are 
retained from the previous example. As is seen, the tempo- 
ral evolution of M strongly depends on the radial distance 
r. At r < 10 4 AU, the mass accretion rate has a well-defined 
maximum at t w 0.21 Myr, when the central gas density has 
exceeded 10 10 cm -3 (the central stellar core formation). Af- 
ter stellar core formation, M drops as ~ t _1 . At r > 10 4 AU, 
the temporal evolution of M does not show a well-defined 
maximum. 



